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Collective response of DNA to THz electric fields is studied in a simple pair bond model. We 
confirm, with some caveats, a previous observation of destabilising DNA breather modes and explore 
the parameter-dependence of these modes. It is shown that breather modes are eliminated under 
reasonable physical conditions and that thermal effects are significant. 
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I. INTRODUCTION 

There is longstanding speculation that nonionising ra- 
diation can damage biological function at the cellular 
level[l]. More specifically, it has also been suggested 
that nonionising radiation of varying frequency causes 
cancer [2|. Since nonionising radiation cannot directly 
disrupt DNA structure, such genotoxic effects must de- 
rive from resonance phenomena driven by external elec- 
tromagnetic radiation^ 

We shall shortly see that the natural frequency of os- 
cillation of DNA base pair separations is approximately 
1 THz, thus THz radiation is of special interest. Further- 
more, it is very likely that this is the unique frequency 
range of relevance to bio-resonance effects in DNA. Inter- 
est in this issue has recently been heightened due to the 
deployment of full body scanners in airports that employ 
millimeter wave (typically 30-300 GHz) technology. 

Motivated by these observations, Alexandrov et al. 
have examined the effects of coupling an electric field 
driving force to a model of dsDNA bond dynamics^. 
The resulting model of damped, driven, coupled non- 
linear oscillators can naturally exhibit exotic collective 
behaviour (for similar earlier conclusions see Ref. [S]). 
For example, the familiar period-doubling approach to 
chaotic dynamics is present. Of more immediate inter- 
est is the discovery of a nonlinear discrete breather mode 
that arises in response to a specific perturbation of the 
system. This mode stores energy for very long times and 
can lead to unbinding effects in dsDNA, with obvious 
implications for the genotoxicity of THz radiation. 

While the results of Alexandrov et al. are compelling, 
it is unclear if the model is sufficiently robust to per- 
mit application to physically realisable situations (such 
as body scanners). In particular, this paper critically ex- 
amines the choice of parameter values, investigates the 
effect of including thermal fluctuations, and examines 
the stability of breather modes in a variety of scenar- 
ios. It will be shown that parameter variation can elim- 
inate breather modes entirely, or make them unrealis- 



tically strong, that thermal noise completely dominates 
the external influences of the system, and that it is ex- 
tremely unlikely that dsDNA denaturing can be induced 
by experimentally accessible THz radiation. 

II. MODEL DEFINITION AND THE AGBUR 
BREATHER 

The model of Alexandrov et al. is based on a model of 
dsDNA pairing dynamics due to Dauxois et al. (the PBD 
model) [6j . The PBD model employs a Morse potential to 
model hydrogen bonding between base pairs (and other 
effects) and an inter-pair stacking potential. Since the 
nucleotide bonding interactions are much weaker than 
those of the phosphate-sugar backbone, the degrees of 
freedom associated with the backbone are neglected. The 
model also ignores degrees of freedom associated with 
the helicoidal structure of dsDNA. The resulting model 
is described by 

mm = -U'M) - W'{y,+^,yi) - W'{y,,y,_^). (1) 

where yi / \/2 is the deviation from equilibrium distance 
of the j'th base pair. The Morse potential is given by 

[/,(y) = A[exp(-a,2/)-l]2 (2) 

and the stacking potential between consecutive base pairs 
is modelled as^ 

W{yi,y^-i) = ]^k{y^-yi_if {1 + pe-xY>{-P{Vi ^ Vi+i)) ■ 

(3) 

In general the parameters can depend on the linked base 
pairs and hence can be labelled etc. Properties of 

this model, including the melting transition, were studied 
inRefs. [ilH]. 

Alexandrov et al. chose to supplement the PBD model 
with periodic driving and frictional terms to model the 
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^ It is worth remembering that biological electric noise generates ^ There is an obvious error in the definition of W in Ref. |4]. This 
internal fields with strengths up to 0.1 V/m[3]. has been repeated in Ref. J^. 
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interactions of dsDNA with an electric field. They state 
that the interactions of the base pairs with an external 
electric field are difficult to model and therefore they as- 
sume a simple harmonic driving force. The additional 
terms are then 



— rrii'yyi + A cos fit. (4) 

Evidently the drag term is not required to model the 
interactions with an electric field, however such a term is 
required to produce collective nonlinear phenomena. Fi- 
nally, Alexandrov et al. assumed a homogeneous poly(A) 
DNA molecule with 64 base pairs. Parameters employed 
were m = 300 amu (this was not specified in Ref. 4] ~ 1 
assume the value given in Ref. |6]), I? = 0.05 eV, a = 4.2 
1/A, k = 0.025 eV/A2, /3 = 0.35 1/A, and p = 2.0. A 
relaxation time typical to water of 7 — 1.0/ps was used. 

As expected, this system displays complicated nonlin- 
ear dynamics. Of particular interest is a breather mode 
found by Alexandrov et al. (which we call the AGBUR 
breather) under a perturbation specified by 

5y,{t) = S{t~to) 0.42cos[^(i-io)] ^(-4 < i < 4)(A) (5) 

at frequency D, — 2.0 THz and with a drive force of 
A = 144 pN. The breather was localised to be approxi- 
mately four base pairs wide and had a maximum ampli- 
tude of approximately 0.3 A. The authors note that fluc- 
tuations like that of Eq. [5] can occur thermally and hence 
transcription and genotoxic effects can be expected. 



A. Breather Characteristics 

In preparation for a detailed examination of these 
claims, we first seek to reproduce the AGBUR breather. 
Solutions were obtained via a microcanonical molecu- 
lar dynamics simulation employing the coupled Runge- 
Kutta (RK4) algorithm. This proved extremely accurate 
(with relative deviations in total energy of order 5 • 10"^ 
over 10 ps) and fast. The Verlet method was also im- 
plemented, yielding results in agreement with RK4, al- 
though less accurate. We follow Ref. [Jl and employ 64 
base pairs with periodic boundary conditions. 

A breather mode was found at = 1.0 THz, some- 
what smaller than the 2.0 THz employed in Ref. |5- Al- 
though it was similar in shape, the maximum amplitude 
of this breather was found to be about 4 A. Note that a 
compression of 0.42 A, such as generated by the pertur- 
bation of Eq. [5j raises the energy of a single bond pair 
by approximately 1.2 eV, which represents an enormous 
insertion of energy. In fact the bond length can then be 
expected to recoil to very large distances, with damping 
supplied by the stacking potential. Thus one anticipates 
that large amplitude breathers, such as found here, are 
to be expected. 



Another breather with a double-lobe structure was 
found at a frequency of fl — 1.5 THz. This novel mode 
is shown in Fig. [l] 
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FIG. 1: (Colour online) A breather mode at f2 = 1.5 THz. 
The compressive perturbation was applied at 2.0 ps. The 
colour scale represents y (A). 

One expects that the particular form of the perturba- 
tion is not important for the formation of breathers. This 
has been confirmed by using a compression of the form 

Sy,{t) = S{t-to)Ype{l<i<n). (6) 

For Yp — —0.42 A, it was found that all perturbations 
with n > 4 generated breathers (with the curious ex- 
ception of n = 8). The double- lobe breather was also 
obtained at = 1.5 THz with this perturbation. At 
n — 5 (and A = 144 pN, = 1 THz) one requires a 
compression of greater than 0.3 A to achieve a breather 
mode. 

The parameters of Alexandrov et al. are not the same 
as those of Dauxois et al.; in particular the value of p was 
changed from 0.5 to 2.0(5]. A run with n = 5, Fp = -0.42 
A, and the PBD parameters^ reveals that the breather 
spreads with time, until the entire DNA molecule melts 
after approximately 140 ps (the same happens with the 
perturbation of Eq. [5| . This is our first indication that 
breather dynamics are subtle and that model results can 
depend crucially on parameters. 

It should be noted that the parameters of Refs. [H E] 
are not universally employed. For example, Barbi et al. 
have developed a similar model that couples base pair 
bond extension to helical twist f 10 . They take D — 0.15 
eV, a = 6.3 1/A, /? = 0.5 1/A, and kp = 0.65 eV/A^. 



3 m = 300 amu, D = 0.04 eV, a = 4.4.5 1/A, A: = 0.02 eV/A, 
/3 = 0.35 1/A, and p = 0.5. 
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We implement this by assuming p — 2.0 and setting 
k — 0.325 eV/A^. Notice that these parameters lead to a 
considerably stiffer collection of nonlinear oscillators. In- 
deed, running with the previous drive parameters {A, f2, 
Yp^ n) reveals that the breather damps out within tens of 
ps. This remains true for all drive parameters that were 
tested. The results make it clear that the relatively large 
stacking interaction disperses the putative breather (see 
Fig. §. 




10 20 30 40 50 60 70 BO 90 100 



t (ps) 



2 

p{AU)= pMB{E)dE ^ —V^e- (7) 
Jau 

where pmb is the Maxwell-Boltzmann energy distribu- 
tion and X — AU/{kBT). The asymptotic form of the 
error function has been used to obtain this result. 

A 0.42 A compression at room temperature yields x — 
46.8 and a probability p ss 10~^°. We now assume that 
n pair bonds must be compressed, 10® base pairs per 
dsDNA, 10^ dsDNA per cell, 10" skin cells^, and a solute 
collision rate of 1.0/ps, to obtain the estimate 

V = 10-^°("~2) 

where V is the probability of obtaining one breather fluc- 
tuation in n base pairs per person per year. We previ- 
ously established that n > A and hence conclude that 
such an occurrence is essentially impossible. 

Motivated by this result and the observation that the 
fluctuation-dissipation theorem links friction with ther- 
mal noise, we have explored the properties of the model 
with the addition of Langevin thermal forcing. Thus the 
term 



FIG. 2: (Colour online) System response to a compressive 
perturbation at 2.0 ps with Barbi parameters. f2 = 3.1 THz. 
The colour scale represents y (A). 



Finally, the effects of allowing two base pair types are 
examined. We model this by setting D = 0.05 eV for 
AT pairs and D = 0.075 eV for GC pairs[n]; all other 
parameters are left at their AGBUR values. We find 
that alternating base pairs or a random configuration of 
base pairs destabilises the putative breather after approx- 
imately 20 ps, again illustrating the fragility the breather 
mode with respect to parameter variation. 



III. THERMAL EFFECTS 

To this point, the genesis of breather modes has re- 
lied on the imposition of a perturbative shock (Eqs. [5] 
[6]) that insert substantial energy into the system. How 
reasonable are these shocks? Presumably they must be 
generated by noise within the cell nucleus. This can be 
due to a variety of biological processes such as cell mem- 
brane activity or by simple thermal fluctuations. 

Here we focus on thermal noise and ask the question: 
how likely is it to perturb a system by Syi? Restricting 
attention to a single base pair, a compression of 0.4 A 
corresponds to an insertion of AU = 1.2 eV of potential 
energy to the system. The probability of such a fluctua- 
tion is 



(9) 



has been added to the right hand side of Eq. [T] We 
assume memory-free noise. The fluctuation-dissipation 
theorem then implies 



{v^it)r],it')) = 2mA3lkBT5{t - t'). (10) 

The molecular dynamics algorithm steps in temporal 
units At and therefore we employ the average noise over 
a time interval (t„, t„ -I- At): 



At 



r]i{t)dt. 



Hence 



,_2\ 2mijkBT 



(11) 



(12) 



Thus average noise forces are chosen from a Gaussian 
distribution 



Pthin) = 



~r/(2(f;")) 



(13) 



* THz radiation is heavily attenuated and only penetrates 1-2 mm 
into the body [13] ■ 
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To test the effect of tliermal fluctuations we revert to 
the homogeneous system with AGBUR parameters and 
driving forces and eliminate the perturbative shocks of 
the previous section. System response was computed 
at a variety of temperatures. At low temperature one 
sees a nonlinear mode with period of approximately 8 
ps that is created by the coupled damped driven oscilla- 
tors. Increasing temperature leads to rapidly increasing 
bond length excursions. This can be quantified by plot- 
ting the distribution of y versus temperature, as in Fig. 
[3j One observes that the distribution is largely invariant 
for ksT < 0.003 eV, fcsT = 0.004 is a transition tem- 
perature, and temperatures greater than 0.005 eV seem 
to yield an invariant distribution for y ^ 1/2 A. Thus it 
appears the dsDNA with AGBUR driving melts at ap- 
proximately 0.004 eV (46K). 

If one were to take this result seriously it would im- 
ply complete chromosomal denaturation in all skin cells 
in the presence of THz radiation. But the previous sec- 
tion warns of large parameter sensitivity in this system 
and one must not arrive at conclusions too hastily. In 
fact using a larger solute relaxation time 7 = 2.0/ps or 
increasing the Morse potential strength to I? = 0.12 eV 
stabilises the system at room temperature. 

It should be noted that there is a subtlety concern- 
ing the thermal properties of the PBD model. The form 
of the interaction implies that the equilibrium configu- 
ration is two widely separated strands. This is reflected 
in the partition function, which necessarily diverges. Of 
course this situation is never realised experimentally be- 
cause the system is embedded in a complex environment 
of approximately /^m extent. In our case, the issue does 
not arise because the molecular dynamics simulation is 
microcanonical and because it only need be run for hun- 
dreds of picoseconds, too short a time scale to probe the 
asymptotic dynamics of the system. More details con- 
cerning the thermal properties of the PBD model can 
be found in Refs. [Hill]- I simply add the observation 
that Barbi parameters yields tiny fluctuations; in fact the 
mean bond extension from equilibrium is (y) ~ 0.03 A 
at 350K, which is inconsistent with the experimentally 
observed melting transition. 



IV. REALISTIC DRIVING 

The AGBUR model of driven dsDNA leads to unrea- 
sonable results in the presence of thermal fluctuations. 
We have already seen how parameter variation can al- 
leviate this problem. But another explanation is pos- 
sible, namely the physical assumptions underlying the 
model could be inaccurate. We follow this by idea by fo- 
cussing on the drive term of Eq. |4] Indeed, this term is 
of immediate concern, since the authors of Ref. [1] state, 
"One complication is that the specific physical nature of 
the interactions between DNA and the THz electromag- 
netic field is not known in detail. ... We will here sim- 
ply augment the PBD [model] to include a drive in the 
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FIG. 3: Distribution curves for y for various ksT with AG- 
BUR driving. 



THz frequency range." In fact, I find that the AGBUR 
breather (at zero temperature, with the AGBUR pertur- 
bative shock) requires A > 140 pN to occur. What is a 
reasonable physical value for the driving force? 

Legal limits on THz radiation power densities range 
from 5-10 mW/cm^|13|. The mean magnitude of the 
Poynting vector {S) = E^/{2fioc) relates^ this to the 
maximum electric field strength, Em- Employing the up- 
per limit gives 



E,n < 30V/m. 



(14) 



Assuming that a nucleotide is singly charged then yields 
a maximum force 



A<4-10~^^N 



(15) 



far smaller than that assumed by AGBUR. In fact the 
driving force is coupled to the base pair displacement and 
therefore depends on the nonuniformity of the field over 
the range (y) . Assuming an incident plane wave (with 
wavevector k) reduces the strength of the drive force by 
a factor of (fc • y) ~ 10~^ hence 



A < 4 • lO^^'^N 



(16) 



But this assumes that the pair bonds all lie in an op- 
timal direction (the force must be along y but is propor- 
tional to k ■ y, thus the pair bond must lie in a plane 
defined by the wave propagation direction and the direc- 
tion of electric field oscillation). In reality, DNA is em- 
bedded in a heavily hierarchical structure, ranging from 



The permeabilities of air and water are essentially equal to that 
of the vacuum. 
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the DNA molecule itself to the chromosome. This effec- 
tively randomises the pair bond direction with respect to 
any external field. Thus one should compute A after av- 
eraging over bond directions. The result is proportional 
to contractions of tensors like Eikj . . .kg which is zero 
to all orders. Thus the driving force relies on remnant 
order in chromosomal structure and A must be much 
smaller than 10~^* N. Finally, approximately one half 
of the incident radiation is reflected [13], cell membranes 
and cytoplasm are extremely efficient at screening elec- 
tric fields, even in the THz regime, and the electric charge 
may be mobile [T5]. All of these effects reduce the cou- 
pling further. One must conclude that the electric field 
driving force is many orders of magnitude smaller than 
that required to generate breather modes. 

V. DISCUSSION AND CONCLUSIONS 

The PBD model of base pair dynamics is sufficiently 
rich that interesting collective behaviour can exhibited. 
Under assumptions concerning drag and drive forcing, 
breather modes can be generated at certain resonant fre- 
quencies. Thus, although this work disagrees on the de- 
tails, it agrees with the main conclusions of Ref. The 
stability of breathers under parameter variation has been 
addressed here. We have seen that changing p from 2.0 to 
0.5 or including thermal noise are sufficient to dissociate 
dsDNA under AGBUR driving. Alternatively, employing 
the stifFer Barbi parameters or allowing for a mixture of 
AT and GC base pairs seems to disallow breather forma- 



tion. 

All of these conclusions are based on drive frequencies 
near the resonant frequency of the system, a drag term, 
and a driving term with a magnitude of approximately 
100 pN. However, it has been argued that the magni- 
tude of the driving term is much smaller than this. The 
physical reason is that the source power is rather weak, 
and DNA is heavily screened from external influences by 
the cell membrane, the cytoplasm, and the nucleoplasm. 
The coupling to electric fields is further reduced by the ef- 
fectively random orientation of a base pair displacement 
vector. The field strength necessary (estimated gener- 
ously) to generate breather modes is approximately 10^ 
V/m, which is much greater than the dielectric break- 
down threshold of air (~ 10^ V/m). Thus it appears that 
the analysis of Refs. [H [7] is not relevant to physically 
realisable situations. 

Although strong THz radiation is artificial, DNA has 
evolved in a noisy electrical and thermal environment, 
and it might be expected that the molecule and the pro- 
cesses in which it takes part will be stable with respect 
to external nonionising radiation. Similarly, one would 
expect that all molecular level biological processes are 
immune to low intensity nonionising radiation; although, 
of course, this speculation needs to be confirmed with 
rigorous experiment. 
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